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Nuclear reaction rates in plasmas depend on the overlap (contact) probability of the reacting 
ions. Path integral Monte Carlo (PIMC) calculations are used here to determine these contact 
probabilities, g{0), for the one component plasma (OCP) with emphasis on many-body quantum 
effects which can lead to order of magnitude changes. An intuitive explanation for these effects 
is presented. The small r behavior of g{r) for quantum systems and the relation to free energies 
is then derived and compared to the path integral results. Going beyond the uniform background 
approximation, electron screening effects and the limits of the "constant energy shift" approximation 
are discussed. Thermodynamic properties for the quantum OCP are analyzed in a final section. 

I. INTRODUCTION 

Due to its importance in calculating nuclear reaction rates, the small r behavior of the radial distribution function, 
g{r), in ionized systems has been a subject of major interest since the original work of Salpeter^. For the uniform 
background model of a classical plasma, simulation studies have produced a quantitative picture of the role of many- 
body effects^. Changes due to quantum effects are less thoroughly understood^*^. 

This article uses path integral Monte Carlo to compute g{r) and in particular 17(0) for dense plasma models. Section 
II presents extensive results for the quantum one component plasma model of ions in a charge neutralizing background 
and gives an intuitive explanation for the trends observed. For classical systems g{0) has a simple connection to a 
free energy difference. The relation is more complicated in the quantum case and section III along with the appendix 
discusses the small r behavior of g{r) in general. Section IV goes beyond the uniform background approximation to 
consider the effects of electron screening^. The thermodynamic functions from the simulations are given in section V. 

The ensuing question of how this information on contact probabilities is ultimately used in a reaction rate calcula- 
tion^ has not been rigorously answered in the literature. The common expression for the equilibrium nuclear reaction 
rate is obtained by multiplying the experimentally determined cross section for a particular reaction, usually as a 
function of incident momentum, by a particle flux and averaging over the Maxwellian distribution for the relative 
incident momentum. Multiplying by the target density gives the reaction rate per volume between species 1 and 2 at 
temperature fcsT = 1/(3 

n„^„,,,l,ip)['-)[J-f'\~^''''^i'p. (1) 

_ For charged particle reactions, the cross section a{p) is usually written as a "nuclear" cross section or astrophysical 
I . factor times a term representing the Coulomb barrier penetration. This last step has the advantage of permitting cross 
''O ' sections measured at accessible experimental energies to be more reliably extrapolated down to the lower energies 
. relevant to most astrophysical applications. 
O ' Different prescriptions have been advanced for modifying this reaction rate to take into account many-body ion and 
electron screening effects. The most common prescription is to simply multiply the above reaction rate by the relative 
change in the contact probability, e^^"-* defined below. For electron screening effects a constant energy shift (see 
section IV) in the Coulomb barrier penetration term is, however, frequently used. The work of Brown and Sawyer— 
perhaps provides a starting point for a fuller discussion of this key question. 



II. OCP CONTACT PROBABILITIES 



Dense plasma effects on nuclear reaction rates are usually described in terms of the enhancement in the contact 
probability 5(0). g{0) is commonly factored into the two-body term, (?6i„(0), and a term, H{0), representing many- 
body effects, 

g{r) ^ gt,Ur)e"^^K (2) 
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For a classical system gbinir) is just e^^^^^^e /r ^ p^^. ^ quantum system gbin{f) is obtained from the solution of the 
Bloch equation for the density matrix. 

For the repulsive Coulomb potential guni^) has an analytic expressioniS 

guniO) = (47r/3) / — g^/fc _ ^ • (3) 

A plot of (76in(0) is shown in the lower panel of Fig. 1. Unlike the classical gbin{r) which is zero at the origin, the 
quantum (7f,m(0) is finite at the origin and increases with temperature. As expected gunif) converges to the classical 
result for r larger than the de Broglie thermal wavelength, = /2'KmkBT . This is demonstrated in the upper 
panel of Fig. 1. 

In this section, contact probabilities are discussed for the OCP model consisting of a single ion species in a uniform, 
charge neutralizing background. In the classical limit, the nonideal properties of this model depend only on the 
dimensionless combination of temperature and density given by the coupling constant F = Z^e^/fc^Ta, where a is the 
ion sphere radius defined by A-Ka^N/V — 1. F is a measure of the relative importance of potential to kinetic energy. 
Contact probabilities for this classical model have been studied starting from the Monte Carlo simulations of Brush, 
Sahlin and Tellcr'^'^^. Empirically H{0) is dominated by a linear dependence on F so it is convenient to define an 
enhancement factor, h{Q) = H{Q)/T. 

When quantum effects are included both density and temperature must be specified and it is necessary to introduce 
a quantum parameter rj = F/r^, with = a/ao, where ao is the Bohr radius for the ions, r] rewritten as ry = 27rA^/a^, 
is seen to be proportional to the squared ratio of the dc Broglie thermal wavelength to the ion sphere radius and thus 
provides an appropriate gauge for quantum effects. 

The many-body g{r) is computed here by averaging over the density matrix e^^^ using Path integral Monte Carlo 
based on the identity 



e 



-f3H ^ ^-I3H/M 



M 



(4) 



where M is an arbitrary integer. Insertion of complete sets of states between the M factors on the right hand side of 
this equation leads to the usual path integral formulation of the density matrix, written here in real space, 

(R|e-''^|R') = p(R, R';I3)= J ...J p(R, Ri; r) . . . p(Rm-i, R'; r)dRi . . . dR^-i (5) 



3 




FIG. 2: Pair correlation functions g{r) and the corresponding many-body enhancement h{r) are shown as a function of the 
couphng parameter, F, for quantum parameter, rj = 1. r is in units of the nuclear Bohr radius. 



with T = P/M . Each of the M steps m the path now has a high temperature density matrix p{R,k, R-fe+i; t) associated 
with it. The integrals are evaluated by Monte Carlo methods. First applied to realistic systems in reference 12, the 
implementation details may be found in recent reviewsi^. In the results presented here the high temperature density 
matrix was taken as a product of exact pair density matrices. Typical M values of 10 to 400, depending on /3, gave a 
discretization error well below the Monte Carlo statistical uncertainties. 

Typical examples of g{r) and h{r) = H(r)/T are shown in Fig. 2. Since g(0) is now finite the h(r) curves 
are considerably easier to extrapolate to the origin than for the classical OCP. A detailed discussion of the formal 
properties of g{r) at small r is given in section III and the appendix. ^(0), gbin{0) and h{0) are given in table I for a 
range of F and rj and are displayed in Fig. 3. 
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TABLE I; Summary of OCP contact probabilities. The coupling parameter F = P/vs and quantum parameter r] — T/rs, or 
equivalently /3 and r^, are listed in columns one through four. The many-body enhancement factor, h{0), is defined from the 
contact probabifity, g{0) = (/6m(0)e^'''-' = g-P(0)+rh(o) (columns five through seven). 
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A clear reduction in /i(0) from the classical value (solid circles), which becomes more important as rj increases, is 
seen in Fig. O This yields reaction rates orders of magnitude smaller than would be predicted by the classical value 
for h{0) at large F. 

Fig-fflprovides the basis for an intuitive understanding of this reduction due to quantum effects. P{r) — — higbini'r) 
can be roughly viewed as proportional to an "effective" quantum pair potential. Increasing rj means that the near 
neighbors of the reacting pair, located approximately one ion sphere radius away, are within a de Broglie thermal 
wavelength. Their effective quantum pair potential is then much less than the Coulomb potential as seen from the 
top panel of Fig. ^ This reduced repulsion lessens the many-body enhancement for the quantum system compared to 
that of a classical system at the same F. Reduction from the classical Coulomb value for h{0) is also seen for classical 
screened Coulomb systemsi^ and for quantum screened Coulomb systems (section IV). The physical explanation is 
again the reduced effective repulsion between the reacting pair and its surrounding neighbors due either to screening 
or quantum effects or both. This is a common pattern. A softer effective potential, either from quantum effects or 
screening enhances the two body contribution to g{0) but reduces the many-body contribution. 

Although /i(0) provides a compact way of presenting the data, the full .g(0), shown in Fig. Ogives a more intuitive, 
physical picture. Fig. 01 shows the transition from the thermonuclear, temperature dependent regime, at low F (high 
temperature) to the pycnonuclear, density dependent ground state regime at large F. The low temperature limit seems 
to be reached when the de Broglie thermal wavelength is of the order of 1/3 the ion sphere radius. The principle 
feature of this graph is the smooth behavior of ^(O). g{0) is seen to be a increasing function of the temperature 
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FIG. 4: Contact probabilities, g{0), versus F for indicated rs showing the transition from thermonuclear (strong temperature 
dependence, low F) to pycnonuclear (temperature independent) regimes. 



and the density. For this model there are no peculiar combinations of density and temperature where the contact 
probabihty, and by implication the reaction rate, has an unexpected local maximum. 
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FIG. 5: Intuitive illustration of the relationship between g{0) and the "mixture" free energy for a quantum system and its 
classical limit. The left panel shows a typical set of discretized (10 steps) paths which would enter the calculation of 17(0) for 
a six particle system. The nodes for the paths of the two "reacting" particles (in center) are shown as open squares. Nodes 
of the surrounding particle paths are shown as filled circles. Because of the delta function in the definition of g{0) the two 
paths for the reacting pair overlap at imaginary times t — and l3. The right panel shows typical paths for the density matrix, 
integrated to get the free energy, of a five particle mixture where the particle in center (nodes shown as open squares) has a 
mass and charge equal to the combined mass and charge of the reacting particles in the left panel. In the classical limit, when 
the de Broglie thermal wavelength goes to zero, the paths for the reacting pair (left panel) or mixture particle (right panel) 
become straight vertical lines. The potential energy at the reacting pair or mixture particle due to the four neighbors is then 
the same for both systems implying that g{0) can be derived from the free energy difference. 



In a classical system the coefficient of the term is proportional to mean squared force on the two fused particles. 
(7(0) is related to the difference between the free energy of a mixture consisting of the fused pair and the N-2 other 
particles and the free energy of the original systemic. For a quantum system the expansion coefficients are more 
complicated. 

Since the delta function in the definition of g{0) causes the two "reacting" particles to overlap it might be thought 
that the relation to the hypothetical mixture free energy which holds in the classical case would be sufficient in general. 
For those comfortable with the path integral ideas used in this paper the cartoon in Fig. 5 may give an intuitive 
understanding of why this is not true. Others are relegated to the detailed derivation in the appendix. 

The two panels of this cartoon show the paths (discretized here into 10 segments) that contribute to g{0) for a 
six particle system (left panel) and the paths contributing to the density matrix of a mixture of five particles (right 
panel) where the two "reacting" particles have been fused into one particle with the combined charge and mass of 
the reacting particles. In the figure the nodes in the discretized paths of the reacting (or fused) pair are shown as 
open squares. The obvious difference in the paths is that in the ^(O) case the delta function in the definition causes 
the paths of the reacting particles to overlap at imaginary times t = and t = (3 but not at other times. In the 
mixture, by contrast, there is only one path for the fused particle. Since this fused particle has the combined mass of 
the reacting pair it would be a "more classical" particle and typically have a less fiuctuating path. 

The scale for the deviation of these paths from straight vertical lines is set by the de Broglie thermal wavelength. 
The different contribution from the two cases (left and right panel) comes from the potential energy around the paths 
of the reacting (or fused) pair due its neighbors. In the classical limit the de Broglie thermal wavelength goes to zero 
and the paths reduce to lines. The potential energy of the other particles at the reacting pair (left panel) is then the 



III. g{r) AT SMALL r AND RELATION TO FREE ENERGIES 



The radial distribution function may be expanded as 
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FIG. 6: Debye-Hiickel approximation for /i(0) (dashed line) as low F limit for the classical, rj — 0, system. PIMC (solid 
diamonds) and semiclassical values (dot-dashed line), Eq. 10, for h{0) at rj — 0.1 are also shown. Hypernetted chain integral 
equation results (crosses) are from Ref. 18 

same as for the fused particle in the mixture case (right panel). Hence the correspondence in the classical limit. This 
also indicates that the classical relation between g{0) and the free energy difference holds when the reacting particles, 
but not necessarily the other neighbors, are treated classically. 

Details for the g(r) expansion are given in the appendix. It follows the original work of Jancovicii^ and Alastuey 
and Jancoviciii. We differ in explicitly separating the gtin (0) term, which is given numerically by Eq. |31 and also 
the free energy difference, approximations for which may be available from other theories. This involves a slight 
rearrangement of terms. 

To lowest order in the quantum parameter 77 and r. 



H{r) = -/3[F(1, N-2)- F(0, N)] - — (r^) 



(7) 



where the first argument in the interaction free energies, F, denotes the number of combined mass, combined charge 
ions, ^r^) is calculated from the relative pair density matrix 



(r2)(r,/3)^- ds 
In the classical limit (r^) ~ and the usual result 



J p(r,r';/3- s)r'2p(r',r;s) 



p(r,r;/3) 



(8) 



H{r) = -/3[F(1, N~2)- F(0, N)] - —r^ 



(9) 



is regained. Using the Debye-Hiickel free energiesi^ in Eq.EI gives the classical weak coupling limit, /i(0) = VST. As 
suggested in Fig. 6 the classical values will approach this for F < 0.1. Although simulation results for the classical 
/i(0) at F < 1 have not been published the approach to the Debye-Hiickel limit is clearly seen in calculations based 
on the HNC approximationi^i. 

More generally Eq. |S1 may be written as 



(r^) (r,/3) = +a(r,/3) 



(10) 



where the defined function a(r, multiplying the free particle result, tends to 1 as /3 ^ 0. a(0,/?) is tabulated in 
Table HI of the appendix. Using Eq. 1101 together with the first order Wigner-Kirkwood quantum correction for the 
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FIG. 7: gbin{0) for the screened Coulomb potential ai k — 100 (open circles), 10 (open squares), and 1 (open diamonds) 
compared to pure Coulomb result (thick solid line) and "energy shift" approximation g5i„(0) = e^''gbin''{0)- 

free energySS, Eq. [71 gives the lowest order quantum correction to the enhancement factor 

hmr,7j] - Mo)[r,o] « -|[a(o,rVr?) - 1/2] . (11) 

The result is shown for 7^ = 0.1 in Fig. 6 and quaUtatively reproduces the PIMC results showing significant reduction 
from the classical h{0). 

IV. ELECTRON SCREENING EFFECTS 

The OCP model discussed above is justified when the electrons are sufficiently degenerate to not respond to the 
ionic potential. Small deviations from this may be treated by linear response leading to an effective ion-ion potential. 
The most common such potential is the screened Coulomb potential 



y(r) = 



with screening length l//t, which we consider in this section. 

From the Feynman-Kac formula applied to the pair density matrix for r = 0, 



Pbm(0,0;/3) 



BMP 



(12) 



(13) 



where the angular brackets denote an average over all Brownian motion paths beginning and ending, after a time f3, 
at the origin. If the screened Coulomb potential is approximated as 
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for r < 1/k then Eg. ll3lbecomes 
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FIG. 8: (color online) h{0) for the classical—, rj = 0, and quantum, 77 = 1, Yukawa system at F = 40 as a function of the 
inverse screening length « times ion sphere radius a. k — corresponds to the unscreened Coulomb potential results in Fig. 3. 



or 5ftj„(0) = e^''g^°"'°'"''(0). The first factor may be interpreted as an energy shift AE = — k. This "constant energy 
shift" approximation fails when the de Broglie thermal wavelength, which gives a scale for the extent of the Brownian 
motion paths, is larger than the screening length, so the paths sample regions where Eg. 1141 does not apply. 
Fig. 7 shows gtiniO) for several k values. These were computed here by rewriting the Feynman-Kac formula as 

^'-£'n m ^(e-fo(ylr(s)]-l/ris))ds\ (^g) 

(0,0;/3) \ /cMP ^ ' 

where now the difference between the screened and unscreened Coulomb potential is integrated over paths distributed 
according to the unscreened Coulomb potential density matrix (CMP). 

The dashed lines show the corresponding constant energy shift approximation. This constant energy shift approxi- 
mation is valid in most astrophysical applications^ but, as seen in Fig. 7, it can dramatically overestimate gbini^) 
low temperatures giving very misleading, excessive reaction rateaSi. 

Turning now to many-body effects, screening has been shown to reduce h{0) in the classical OCP— . Its effect in 
the quantum OCP, shown in Fig. 8, is similar. The reduced repulsion from surrounding ions due to screening again 
reduces the enhancement effect. 



V. QUANTUM OCP THERMODYNAMIC PROPERTIES 

The PIMC calculations of the density matrix used to obtain g{r) also yield the kinetic and potential energies. 
Results for the computations used in this paper are given in Table II. 
The PIMC and harmonic approximation values for the kinetic energy 
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where ujjCk) are the vibrational frequencies of the BCC Wigner lattice24 are compared in the table. Even though all 
simulations tabulated were in the fluid phase the agreement at higher F is quite good. This is not entirely surprising 
as similar agreement with other properties has been often noted for the purely classical system due to the long range 
nature of the interaction. The harmonic approximation also gives the correct lowest order quantum correction. 
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FIG. 9: (color online) Comparison of the quantum kinetic, K — Kc (open symbols and dashed lines), and potential energy, 
V — Vc (solid symbols and lines) per particle (data from Table II). 




FIG. 10: (color online) Comparison of the quantum energy per particle, E — Ec for F = 1 (solid circles), F = 2 (solid squares) 
F = 5, (solid diamonds), F = 10 (crosses) and F = 40 (solid triangles) with the second order Wigner-Kirkwood {h expansion) 
correction (corresponding solid or dashed lines). 



using the Kohn sum rule. As expected for fixed F this agreement worsens as rj increases. 

The quantum corrections to the kinetic and potential energy as functions of F and r] are plotted in figure 9. As 
expected the magnitude of the quantum corrections for both quantities increases with i]. At high F the quantum 
kinetic and potential energies are seen to converge. Again, if the harmonic approximation where theses quantum 
corrections are equal, gives a good description of the thermodynamics even in the liquid state at high F, then this 
convergence is understandable. The slower convergence as rj increases for a fixed F is also understandable since 
anharmonic corrections would be larger here. 

The quantum contributions to the total energy are compared to the expression for the second order Wigner Kirkwood 
correction as given by Hansen and Viellefosse^^ in Table 11 and in Fig. 10. Similar limits on the convergence of this 
second order approximation were previously seen by Jones and CeperleyM. (See their figure 3.) 
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TABLE II: Excess quantum kinetic {K — Kc), "excess" potential (V — Vc), and internal (E) energy per particle from PIMC with 
54 distinguishable particles and M time slices. Kc ~ 3/2/3 and Vc is the potential energy for the classical OCP with 54 particles 
in the periodic ceU, taken as /3K = {-0.24497,-0.57994,-1.32795,-3.7621,-7.9997,-34.2373} for T = {0.5,1,2,5,10,40} 
and as given by the fit^'^ /3Vc = -0.899375r + 0.569333r^/^ - 0.224470 - 0.017875/r^/^ for larger T values. The subscript h 
denotes values from the harmonic approximation for the BCC Wigner lattice. The first and second order Wigner Kirkwood 
corrections to the energy from reference 25 are given in the last column 
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0.01766(8) 


0.00691(8) 


0.00713 


0.25 


5 


100 


10 


80 


0.00138(8) 


0.0015 


-0.036220(4) 


0.001401(4) 


-0.01984(8) 


0.00278(8) 


0.002819 


0.25 


10 


400 


40 


40 


0.000791(7) 


0.000727 


-0.0193250(6) 


0.0006743(6) 


-0.014784(7) 


0.001465(7) 


0.001348 


0.25 


40 


6400 


160 


80 


0.000146(2) 


0.000156 


-0.0052058(1) 


0.0001438(1) 


-0.004826(2) 


0.000290(2) 


0.000241 


0.25 


100 


40000 


400 


200 


0.000047(2) 


0.000051 


-0.00213982(9) 


0.0000472(1) 


-0.002056(2) 


0.000094(2) 


0.000022 


0.50 


0.5 


0.5 


1 


50 


0.015(8) 


0.016 


-0.4389(1) 


0.05099(1) 


2.576(8) 


0.066(8) 


0.06798 


0.50 


1 


2.0 


2 


50 


0.014(2) 


0.0306 


-0.26155(4) 


0.02842(4) 


0.503(2) 


0.0424(2) 


0.04797 


0.50 


2 


8.0 


4 


40 


0.0109(4) 


0.0151 


-0.15176(1) 


0.01423(1) 


0.0466(4) 


0.0251(4) 


0.02582 


0.50 


5 


50 


10 


50 


0.00491(9) 


0.00582 


-0.069885(4) 


0.005356(4) 


-0.03497(9) 


0.01027(9) 


0.01005 


0.50 


10 


200 


20 


50 


0.00231(2) 


0.00274 


-0.037492(2) 


0.002506(2) 


-0.02768(2) 


0.00482(2) 


0.00453 


0.50 


40 


3200 


80 


80 


0.000593(5) 


0.000541 


-0.0102038(6) 


0.0004954(6) 


-0.009142(5) 


0.001088(5) 


0.000364 


0.50 


100 


20000 


200 


100 


0.000145(2) 


0.000166 


-0.0042185(2) 


0.0001555(2) 


-0.003998(2) 


0.000301(2) 


-0.000452 


1.0 


0.5 


0.25 


0.5 


25 


0.066(9) 


0.245 


-0.8154(2) 


0.1644(2) 


5.251(9) 


0.230(9) 


0.04381 


1.0 


1 


1.0 


1 


10 


0.041(1) 


0.121 


-0.48277(6) 


0.09722(6) 


1.059(1) 


0.138(1) 


0.1338 


1.0 


2 


4.0 


2 


40 


0.027(2) 


0.059 


-0.28043(5) 


0.05156(5) 


0.121(2) 


0.079(2) 


0.0816 


1.0 


5 


25 


5 


25 


0.0141(2) 


0.0219 


-0.13102(2) 


0.01947(2) 


-0.0569(2) 


0.0336(2) 


0.03044 


1.0 


10 


100 


10 


80 


0.00723(9) 


0.0100 


-0.071170(6) 


0.008827(4) 


-0.04894(9) 


0.01606(9) 


0.0112 


1.0 


40 


1600 


40 


128 


0.00156(2) 


0.00178 


-0.019772(2) 


0.0016270(7) 


-0.01727(2) 


0.00319(3) 


-0.00334 


1.0 


100 


10000 


100 


100 


0.000495(2) 


0.00052 


-0.008258(2) 


0.0004905(7) 


-0.007606(2) 


0.000986(3) 


-0.00611 


1.0 


200 


40000 


200 


400 


0.00021(2) 


0.000198 


-0.004248(3) 


0.000172(3) 


-0.00399(2) 


0.00038(2) 


-0.00697 


1.0 


400 


160000 


400 


400 


0.000068(4) 


0.000074 


-0.0021575(9) 


0.0000662(9) 


-0.002080(4) 


0.000134(4) 


-0.00735 


1.0 


600 


360000 


600 


450 


0.000047(6) 


0.000041 


-0.001453(2) 


0.0000334(2) 


-0.001402(7) 


0.000080(7) 


-0.00746 


2.0 


0.5 


0.125 


0.25 


50 


0.14(3) 


0.97 


-1.4797(4) 


0.4799(4) 


10.66(3) 


0.62(3) 


-1.6495 


2.0 


1 


0.5 


0.5 


50 


0.08(2) 


0.471 


-0.8625(3) 


0.2974(3) 


2.22(2) 


0.38(2) 


0.07017 


2.0 


2 


2.0 


1 


50 


0.074(4) 


0.224 


-0.4993(1) 


0.1647(1) 


0.325(4) 


0.239(4) 


0.1526 


2.0 


5 


12.5 


2.5 


50 


0.035(1) 


0.080 


-0.23622(5) 


0.06475(5) 


-0.081(1) 


0.010(1) 


0.04349 


2.0 


10 


50 


5 


50 


0.0184(3) 


0.0345 


-0.13097(2) 


0.02902(2) 


-0.0826(3) 


0.0474(3) 


-0.01006 


2.0 


40 


800 


20 


80 


0.00387(5) 


0.00565 


-0.037748(4) 


0.005049(4) 


-0.03201(5) 


0.00892(5) 


-0.0517 


2.0 


100 


5000 


50 


100 


0.00145(2) 


0.00158 


-0.016003(2) 


0.001493(2) 


-0.01425(1) 


0.00294(1) 


-0.0589 



VI. CONCLUSIONS 

In conclusion, quantum effects have been shown to significantly reduce the many-body enhancement factor which 
influences nuclear reaction rates in dense plasmas. Electron screening effects produce a further reduction. The contact 
probability and reaction rates based on it increase monotonically with temperature and density. The relation between 
the contact probability and free energy differences is derived and intuitively illustrated. 

APPENDIX A: SMALL r EXPANSION OF h{r) IN THE SEMICLASSICAL LIMIT 

In a classical plasma g{0) can be related to a free energy difference^ and a simple, explicit value given for the 
coefficient of the expansion term-S,. For the quantum plasma no similar, exact relation has been found however 
an expansion in the ratio of the de Broglie thermal wavelength to the ion sphere radius for g{r) at small r can be 
madei^. This expansion, which reduces to the classical result, is reviewed here. The principal results are given below 
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in Eas. rn2l HTfl HTSl and in Table III for a(0;/3). 

Starting from the definition of the radial distribution function for a one component system 

.g(r) = f7(<5(r-ri2)) = 1^ ^^-UH^ ^ ^ ^ 

where particles one and two have been singled out, the Hamiltonian 

N 

j=l k j<k 

with K the kinetic energy operator is rewritten in terms of center of mass and relative coordinates for the " reacting" 
particles, one and two, R = (ri + r2)/2 , ri2 = ri — r2. Doing this, separating out the terms involving only particles 
one and two, and adding and subtracting a term representing the interaction of both particles one and two at the 
center of mass position with the "spectator" particles 3 . . . iV the Hamiltonian can be rewritten as, 

i/ = i/,ei(r)+iJ™„(R,r3...rjv) + AT/(R,r,r3...rAr) . (A3) 

The "relative" Hamiltonian 

Hrei{r) =^ Kr + v{r) (A4) 

commutes with the "mixture" Hamiltonian 

N N N 

iJ„„(R, rg . . . rjv) = i^fl + E ^0 + E E "j'' + 2 E " "^j) 

j=3 j=3 k>j j=3 

which corresponds to a particle of double the mass and charge at the center of mass position, R, and the N — 2 
spectator particles. The coupling term 

AF(R, r, ra . . . r^v) = E'Is H^i - ^j) + «(r2 - rj) - 2v{R - rj)] 

= EL [«(R + § - rj) + i;(R - § - rj) - 2«(R - rj)] (A6) 
= E ^3 VVw(R - rj) : rr/4 + O(r'i) 

is the difference between the interactions of all other particles in the system with particles 1 and 2 at their actual 
positions minus these interactions when particles 1 and 2 are fused at their center of mass position. 
Using Eg. IA3l in the expression for g{r) and taking the trace in real space 

5(r)= I J (R,r,r3,...rN|e-«^-+^— +^^)|R,r,r3,...rN)dRdr3...drN . (A7) 
When averaged over an isotropic system 

VVw(R-rj) = v2i;(R-rj) T . (A8) 

The Laplacian is easily evaluated for the Coulomb system considered here where the interaction, accounting for periodic 
boundary conditions and charge neutrality, is the Ewald potential, "ifEwaidir)- Using W^'^ Ewaid{r) = — 47r5(r) + 47r/r2, 
where, physically, the constant comes from the neutralizing background. 

Ay (R, r, r3 . . . r^.) = ^^^^^^ ^ + 0{r') = ^ + 0{r') . (A9) 

The simplification that the coefficient of the lowest order term in r is constant is unique to Coulomb systems. 
With this expansion and using the fact that Hmix commutes with Hrei and the above lowest order term for AV 

g{r) - I /r|e-/'(^'-+'^'-')|r\/(R,r3,...rN|e-/'«— |R,r3,...rN)dRdr3...drN + 0(<r4 >) 
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where C = Z'^/Aa^. The free energies Fp^re, Z = exp(— /?Fpnre), and Fmix^ 

e-'^^— = y" (R, ra, . . . rN|e-^^""== |R, rg, . . . tn) dRdrs . . . drN (AlO) 

correspond to the original, total Hamiltonian and to Hmix respectively. They are the fully quantum mechanical free 
energies. 

Introducing the binary radial distribution function, from the definition Eg. lAll applied to a two particle system, 

gu^r) = {Ai^f^aimf (r|e-'^^- |r) (All) 

the expansion Eo. IA10| for g(r) becomes 

The term following cancels the remaining "ideal gas term" in the difference between the iV particle Fp^re and 
the iV — 1 particle Fmix- In the following this term is omitted and the free energies refer only to the nonideal terms. 
For the classical case Hj-ei and Cr^ commute and the last ratio in equation above is just e~'-^^ so 

g[r) = g,„(r)e-^(^""-^-"-'^'-^) = g,,„(.)e-'5(^"--^--)-5 ('■/-)^ . (A13) 

This relation, with the correct term, was first derived in reference 16 correcting an earlier factor of two error in 
reference 2. 

As a simple application consider the weak coupling, Debye-Hiickel limit. Using the well known interaction free 
energy for this modcl^^^ 

applied to the "pure" {N^ = N , Zi = Z) and the "mixture" case {Ni ^ N - 2, Zi = Z, N2 = 1, Z2 = 2Z) gives 
h{0) = ln[g{0)/gbin{0)] — VST- This result also comes from expanding the screened potential form, g{r) — e'^^/r 
with K the inverse Debye-Hiickel screening length. The next term in the expansion, linear in r, is seen to be incorrect 
however. 

In the quantum mechanical case the terms Hrei and Cr^ do not commute. The identity 







leads to 



_g-/3ff„, _/ e~('3-^)^"'Cr2e-^^"'ds + 0(r4) (A16) 




j.|g-/3(ff,el+Cr")| 



= 1-I3C <r^ > +0{< >) (A17) 
(r|e-P""'|r) 

where 

(r^) = ds J dr'preiir,r';f3- sy^prei{r',r;s)/prei{r,r;(3) . (A18) 

Unlike the classical case < >^ but, because of the range of prei, will differ from it by an amount proportional to 
the squared de Broglie thermal wavelength. For example, if prei is approximated by only the free particle term then 

< + h^(3/M . (A19) 

When r — > 0, < > and higher order terms are thus nonzero and the simple, classical relation between the screening 
function at r = and the free energy difference no longer applies. What has been generated is a double expansion in 
and the squared ratio of the de Broglie thermal wavelength to the ion sphere radius. 
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TABLE III: a(0,/3) as defined in Eg. 1X221 and IXTsl 



/3 


q(0,/3) 


0.0 


1.0 


0.5 


1.106 


1.0 


1.148 


2.0 


1.206 


4.0 


1.286 


6.0 


1.345 


8.0 


1.394 


10.0 


1.435 


20.0 


1.591 


30.0 


1.703 



The lowest order quantum correction to the free eneigy^ 

I3F = PF,^ass^cal + ^ ^ ^ (A20) 

where U is the total potential energy function, applied to the uniform background model and following the algebra 
leading to equation IA9I gives 

P - = P {F^r^'-' - F^tT""') - . (A21) 

The quantum correction in this term would further increase the enhancement factor contrary to what is found. 

The f3C < > term in Eq. IA17I corrects this. We have not found a simple expression for < > in terms of 
continuum Coulomb wave functions however it is not difRcult to evaluate numerically using the axial symmetry of 
the r' integral (or radial symmetry when r = 0). The result can be expressed as 

{r^){r,l3)^r' + a{r,f3)^ (A22) 

where the function multiplying the free particle result, a(r, /?) — > 1 as /3 — > 0. As /3 increases a slowly increases, 
reflecting the tendency of the repulsive Coulomb potential to emphasize larger radius " paths" compared to the free 
particle limit. 

Adding this term to the free energy change, the lowest order quantum correction in the enhancement factor is 

AhQMiO) - -|[a(0,rV,7) - 1/2] . (A23) 

This now correctly predicts the decrease in h{0) due to quantum effects. a(0,/3) is tabulated in table III. This 
semiclassical expansion was compared with the PIMC results for the case r/ = 0.1 in Fig. El For values of r/ > 0.25 
the lowest order expansion overestimates the quantum effects by almost a factor of 2. 

We thank H. DeWitt for numerous conversations and for several of the classical h{Q) values displayed in figure 3. 
This work was performed under the auspices of the U.S. Department of Energy by University of California LLNL 
contract No. W-7405-Eng-48 in a CREM free environment. 
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